function [TT, TT_deriv_1, TT_deriv_2] = construct_chebyshev_twodim(x_interest_1,x_interest_2,K1,K2,a1,a2,b1,b2)

[T_1, T_deriv_1]=construct_chebyshev(x_interest_1,K1,a1,b1);
[T_2, T_deriv_2]=construct_chebyshev(x_interest_2,K2,a2,b2);

TT=[];
TT_deriv_1=[];
TT_deriv_2=[];

for ii=1:(K2+1)

    T_add=T_1.*(T_2(:,ii)*ones(1,1+K1));
    TT=[TT T_add];    

    T_deriv_add1=T_deriv_1 .* (T_2(:,ii)*ones(1,1+K1));
    TT_deriv_1=[TT_deriv_1 T_deriv_add1];
    
    T_deriv_add2=T_1 .* (T_deriv_2(:,ii)*ones(1,1+K1));
    TT_deriv_2=[TT_deriv_2 T_deriv_add2];
    
end

    %
